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Abstract. We compare velocity structure observed in the Polaris Flare molecular cloud at scales ranging from 
0.015 pc to 20 pc to the velocity structure of a suite of simulations of supersonic hydrodynamic and MHD 
turbulence computed with the ZEUS MHD code. We examine different methods of characterising the structure, 
including a scanning-beam method that provides an objective measurement of Larson's size-linewidth relation, 
structure functions, velocity and velocity difference probability distribution functions (PDFs), and the A- variance 
wavelet transform, and use them to compare models and observations. 
, The A- variance is most sensitive to characteristic scales and scaling laws, but is limited in its application by 

a lack of intensity weighting so that its results are easily dominated by observational noise in maps with large 
empty areas. The scanning-beam size-linewidth relation is more robust with respect to noisy data. Obtaining the 
global velocity scaling behaviour requires that large-scale trends in the maps not be removed but treated as part 
of the turbulent cascade. We compare the true velocity PDF in our models to simulated observations of velocity 
centroids and average line profiles in optically thin lines, and find that the line profiles reflect the true PDF better 
unless the map size is comparable to the total line-of-sight thickness of the cloud. Comparison of line profiles to 
velocity centroid PDFs can thus be used to measure the line-of-sight depth of a cloud. 

The observed density and velocity structure is consistent with supersonic turbulence with a driving scale at or 
above the size of the molecular cloud and dissipative processes below 0.05 pc. Ambipolar diffusion could explain 
the dissipation. Over most of the observed range of scales the velocity structure is that of a shock-dominated 
medium driven from large scale. The velocity PDFs exclude small-scale driving such as that from stellar outflows 
as a dominant process in the observed region. In the models, large-scale driving is the only process that produces 
deviations from a Gaussian PDF shape consistent with observations, almost independent of the strength of driving 
or magnetic field. Strong magnetic fields impose a clear anisotropy on the velocity field, reducing the velocity 
' variance in directions perpendicular to the field. 
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1. Introduction 

Attempts to characterise the physical state of molecular 
clouds by comparison to simulations must rely on sta- 
tistical descriptions of the observations and the simula- 
tions. In the last decade, several techniques have been 
used to characterise the observed radial velocity distri- 
bution in representative molecular clouds, as reviewed by 
Goodman et al. (1998)| and Miesch et al. (1999). 

First attempts to characterise the scaling behaviour 
of the velocity structure started from one of the fa- 



mous "Larson's laws" (Larson 1981), showing a power 
law relation between the size and the linewidth measured 
for a molecular cloud. This relation has been extended 
from inte grated velocities to velo city fluctuations within 
clouds by Miesch fc Bally (1994) , providing similar power 
laws. On the other hand hydrodynamic and magneto- 



hydrodynamic (MHD) simulations have been tradition- 
ally characterised using the probability distribution func- 
tions (PDFs) of velocities and velocity differences (e.g. 



by Lis et al. (1996' 



Ansclmct et al. 1984 ). Bo th approaches have b een unified 
and 



Miesch et al. (1999) , who mea- 



sured the PDFs of line centroid velocities and the scal- 
ing behaviour of PDFs of centroid velocity differences as 
a function of lag for several star-forming clouds. 

Mac Low & Ossenkopf (2000, hereafter paper I) used 
the A-variance, a multi-dimensional wavelet transform 
( Stutzki et al. 1998| ), to characterise both the density and 
the velocity structure of interstellar turbulence simula- 
tions. For the density structure, a direct comparison to 
the analysis of observed clouds provided by Bensch et al. 
(2001a) was possible. In velocity space there is no direct 
observational measure for the A-variance. We need addi- 
tional tools for the quantification of the turbulent velocity 
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structure that can be determined with the same ease for 
observations and simulations, and containing at least as 
much information as the A- variance. 

In this paper we test five different methods on an ob- 
servational data set covering three steps of angular reso- 
lution in the Polaris Flare, a translucent molecular cloud, 
and on a number of gas dynamical and MHD simulations 
of interstellar turbulence. The first two of the methods 
characterise the total velocity distribution: the PDFs of 
the total velocity distribution and of the line ccntroid ve- 
locities. The other three methods characterise the spa- 
tial distribution of velocities: a generalised size-lincwidth 
( "Larson" ) relation, the dependence of the low order mo- 
ments of the centroid velocity difference PDF on lag, and 
the A-variance analysis. 

Comparisons between observations and models have 
been made with a simulation of mildly supersonic, de- 



caying hydrodynamic turbulence ( 


Falgarone et al. 1991, 


Falgarone et al. 1994 




Falgarone et al. 1995, Falgarone 



et al. 1995, Lis 



Joulain et al. 1998 



et 

and 



al. 1996, Lis et al 1998, 
Pcty fc Falgarone 2000[), with 



MHD mo dels of supersonic turbulence neglecting self- 
gravity by Padoan et al. (1998) , Padoan et al. (1999) and 



Padoan et al. (2000) , with decaying and driven self- 
hydrodynamic turbulence (Klessen 2000), 
various ad hoc models of turbulence 
Chappell & Scalo 1999). 
Observations suggest that supersonic, super-Alfvenic 
turbulence is an appropriate physical model (see 



gravitating 
and with 

(e 



Dubinski et al. 1995 



Padoan fc Nordlund 1999 ). Here, we want to test whether 
it can reproduce the observed velocity structure and 
find what we can learn about the physical conditions of 
turbulent molecular clouds from comparison to a large 
set of models computed by Mac Low et al. (1998) and 



Mac Low (1999) 



In Sect. 2 we discuss the Polaris Flare observational 
data used here and the basic data processing applied. In 
Sect. 3 the different tools for the analysis of the velocity 
structure are introduced and applied to the observational 
data. The turbulence simulations used for comparison are 
presented in Sect. 4 and the results of the velocity analysis 
for these models are given in Sect. 5. Sect. 6 concludes with 
a discussion on the physical state of the cloud based on 
our comparisons. 



2.1. Polaris Flare observations 

In Paper I, we compared our analysis of simulations with 
an analysis of the intensity structure of multiscale obser- 
vations of the Polaris Flare performed by Bensch et al. 
(2001a). As no prior analysis of the velocity structure ob- 
served in this data has been done, we present that here as 
one point of comparison to the simulation results. 

For the velocity field analysis we use three of the ob- 
servational data sets studied by Bensch et al. (2001a) in 
the analysis of the intensity structure. The Polaris Flare 
observations consist of a set of nested maps obtained with 
the 1.2 m CfA telescope, the 3 m KOSMA, and the 30 m 
IRAM. The CfA data were taken in 12 CO 1-0 at a spa- 
tial resolution (HPBW) of 8.7' (Heithausen & Thaddeus 
1990); the KOSMA observations used the 12 CO 2-1 transi- 
tion at 2.2' resolution (Bensch et al. 2001a); and IRAM ob- 
servations of the MCLD 123.5+24.9 region in the Polaris 
Flare were taken in the two lower transitions of 12 CO and 
13 CO and in the 1-0 transition of C 18 within the IRAM 
key-project "Small-scale structure of pre-star-forming re- 
gions" (Falgarone et al. 1998). To discuss a consistent set 
of observations for all maps we restrict ourselves to the 
12 CO IRAM data. Because of the higher signal-to-noise 
ratio, we only use the 1-0 lines taken at 0.35' resolu- 
tion. Assuming a distance to the Polaris Flare cloud of 
150 pcQ, the telescope resolutions translate into physical 
resolutions of 0.38 pc, 0.09 pc, and 0.015 pc respectively. 
Altogether these observations provide a data set cover- 
ing more than three decades in linear resolution - from 
0.015 pc to about 50 pc. 

A major problem when combining these data is the dif- 
ferent channel width, noise, and baseline behaviour of the 
different instruments and observational runs. The IRAM 
observations show an rms noise of 0.5 K at 0.05 km/s, 
with some large scale trends in the noise indicating that 
either the spectral or the spatial baseline is not optimal. 
The KOSMA data have an rms of 0.4 K at 0.05 km/s 
throughout the whole map. The CfA data show an rms 
noise of 0.1 K at 0.65 km/s, and show both a variation of 
the absolute noise level throughout the map, and some ar- 
eas with slightly negative intensities, indicating imperfect 
baselines. This opens up some uncertainties when combin- 
ing noise sensitive results from the three data cubes. 



2. Observational data 

Molecular line observations only determine line profiles, 
which give the convolution of the radial velocity compo- 
nent with density along the line of sight. The situation 
becomes even more complicated if one takes optical depth 
effects and spatially varying temperatures and excitation 
levels into account. For the analysis provided here we re- 
strict ourselves to the assumption of constant excitation 
conditions in an optically thin medium, so that the in- 
tegrated line intensity is a direct measure of the column 
density. 



2.2. Line windowing 

The basic problem in the deduction of the velocity struc- 
ture is the finite signal-to-noise ratio in each observed 
map, often combined with an imperfect spectral baseline. 
Due to small slopes and variations of the baselines, and 
a slightly variable noise throughout the spectrum, the de- 
termination of the velocity structure from the spectra is 
sensitive to the exact selection of the spectral window 
around the line considered. The influence of this effect 



|Heithausen fc Thaddeus (1990) 



estimated a distance to 



the Polaris Flare cloud of 240pc whereas Zagury et al. (1999) 
derived 105-125pc. 
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increases from the line centroid velocities to the higher 
moments like the variances or the kurtosis. A detailed dis- 
cussion of these problems in the determination of the cen- 
troid velocity was provided by Miesch & Bally (1994) and 
Miesch k Scalo (1995) . 

We have tested three different methods for the selec- 
tion of that part of the spectrum containing as much infor- 
mation as possible about the line but least influenced by 
baseline uncertainties. First we applied a global window- 
ing technique defining a minimum velocity range covering 
all noticeable emission. We used this global windowing as 
a first rough constraint to the velocity space, including 
all channels where eye inspection might still guess some 
line contribution. A second method that we have exten- 
sively tested was to search for the line contribution in- 
terval by using the first zeros at the flanks of the lines. 
We found, however, that too many lines are either rela- 
tively weak or break up into several components so that 
zeros occur within the line. Even the centroid velocities 
could not be reliably determined from this approach. As 
a third approach, we have used a typical criterion for no- 
ticeable emission, such as emission above a 3 a noise level. 
It turned out, however, that for most values of the signif- 
icance level, we clearly miss part of the information from 
weak but broad lines in thin outer regions that inspection 
by eye would still count as a part of the line. 

Thus we show for all our observational results large er- 
ror bars given by two extremes of this criterion: we count 
cither all contributions within the global spectral window, 
or only channels above the maximum noise level given by 
the largest negative value. As representative intermediate 
values we show in all plots the results using all positive 
contributions above la of the noise. This third level pro- 
vides a reasonable intermediate value and the parameters 
from the data analysis below show only a small variation 
when changing the noise cut level around this value. One 
should however keep in mind that we do not know the 
best treatment of the noise so that this line should be 
considered only as guiding the eye but not as the best 
representative or average value. 

2.3. Large-scale trends 



Miesch fc Bally (1994J extensively discussed the removal 
of large scale trends in the centroid velocity maps to get 
a significant description of the turbulent velocity fluctu- 
ations undisturbed by any large systematic motions. We 
will not follow their method here when dealing with the 
Polaris Flare data. Due to the nested nature of the dif- 
ferent maps, any large scale motion on a smaller map 
is only a velocity fl uctuation on the larger map. As 



Bensch ct al. (2001a) have shown for the intensity maps, 



a smooth transition in the scaling behaviour of the three 
maps is only possible if large scale trends are not removed. 

Furthermore, there is no clear separation between tur- 
bulent and systematic motions. Any assumed separation 
scale is arbitrary as long as there is no physical process 



such as energy injection at that scale. Removal of velocities 
at certain scales might prevent understanding of the un- 
derlying processes. Even the largest systematic motions, 
like Galactic rotation, may be part of the turbulent cas- 
cade if they inject energy into the system. 

It could be justifiable to remove the large scale 
trends for the star-forming regions considered by 



Miesch fc Bally (1994) if the turbulence there were only 
driven by small-scale star-formation activity, and not by 
large-scale motions. However, that is not proven even 
there, and we would certainly miss the main physics by 
applying the same kind of separation for the Polaris Flare 
data, where the turbulence probably is driven by motions 
on the largest scales (see discussion in Paper I). 

3. Statistical descriptions of velocity structure 
and their applicability to observations 

3.1. Size-linewidth relation 

A traditional measure for the spatial velocity distribu- 
tion is the size-linewidth relation for clouds identified by 
Larson (1981) and obtained by many observers since then. 
To measure this, one has to define objects within a map, 
such as molecular clouds or clumps within clouds, and re- 
late the effective linewidths measured for these objects to 
their characteristic sizes. This method has been used to 
study a large variety of clouds, clumps, and cores (e.g. 
paselli fc Myers 1994 ~ 



Myers 1983 



Peng et al. 1998). A 



comprehensive recent overview including a careful esti- 
mate of many possible errors was given by Goodman et 
al. (1998). Most studies obtain power laws 



Av ohs (x W with 7 = 0.2 . . . 0.7 



(1) 



over wide spatial ranges, where R denotes the effec- 
tive radius of the object. A search for Larson-type 
relations in turbulence simulations was performed by 
Vazquez-Semadcni ct al. (1997)[ 

A major problem in the computation of these size- 
linewidth relations is the somewhat arbitrary definition of 
the objects in the observed position- velocity space that are 
considered to give definite values for sizes and linewidths. 
This definition is obvious, though dynamically arbitrary, 
for isolated molecular clouds, but difficult when selecting 
clumps within a cloud, and completely impractical for fil- 
amentary, turbulent cloud structures. 

As an alternative, one can compute a size-linewidth 
relation by measuring the average linewidths within tele- 
scope beams of varying size - effectively the intensity- 
weighted velocity dispersion within a varying radius. 
Practically, the observed map is scanned with a Gaussian 
of varying size and the average linewidth is determined for 
each size. In averaging, each position is weighted by its 
total intensity so that the lower significance of "empty" 
regions is taken into account. This method can be eas- 
ily applied to each data set and the results are directly 
comparable to the traditional size-linewidth relations. 
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Fig. 1. Sizc-lincwidth relation for the Polaris Flare CO 
observations with IRAM (smallest scales), KOSMA, 
and the CfA 1.2m telescope (largest scale, see 
Bensch ct al. (2001a)| for details) . The diamonds show the 
relation when the linewidth for a given scale is integrated 
from the total linewidths at each point. The triangles rep- 
resent the widths when only the dispersions of the cen- 
troids of the lines are measured. The error bars do not 
represent true local errors but the two extreme cases of 



the line windowing as discussed in Sect. 2.2 



We face, however, the problem that the velocity disper- 
sion for each virtual beam depends not only on its size, but 
also on the cloud depth along the line of sight. This depth 
enters by changing the local linewidth at each point. To 
separate the two length scales we have applied the analysis 
both to the total velocity dispersion within the Gaussian, 
and to the dispersion of the line centroids observed at each 
point in the map. In the latter case, we hope to remove 
the influence of the cloud depth. 

In Fig. [l] we show both relations for the Polaris Flare 
observations. The two extreme approaches for the noise 
treatment (taking all emission from the spectral window 
or only emission above the maximum noise level) produce 
relatively large error bars in the plots for both the total 
velocity dispersion and the centroid velocity dispersion. 
The CfA data with their low spectral resolution in par- 
ticular have error bars of up to a factor of two. Within 
the errors, however, we find a unique smooth behaviour 
in both quantities from the smallest to the largest scales. 
The KOSMA data do deviate somewhat from the other 
two sets. This can be partially attributed to the differ- 
ent CO transition observed by that telescope. We have 
checked this effect by computing the same plot for the 
IRAM 12 CO 2-1 data. Although the results were shifted 
relative to the IRAM 1-0 data in the same direction as 
the KOSMA data, they do not line up exactly with the 
KOSMA result. Thus the shift is probably also influenced 
by the different noise behaviour. 

For the size-linewidth relation based on the velocity 
centroids, we find one power law stretching over three or- 
ders of magnitude connecting the three different maps. 



The average velocity variances range from below the ther- 
mal linewidth up to about 1 km s . The common slope is 
given by 7 = 0.50 ± 0.04. However, the data are also con- 
sistent with a reduction of the slope down to 0.24 at the 
largest scales, if the full extent of the error bars is taken 
into account. 

In the size-linewidth relationship integrated from the 
full local linewidths, there is a transition of the slope from 
almost zero at scales below 10' to 0.2 at the full size of 
the flare. The plot shows that the total linewidths are 
dominated by the line-of-sight integration up to the largest 
scales. 

Although the slopes measured with this method are 
very shallow, they do appear to show the change of slope 
interpreted by Goodman ct al. (1998) as a transition to 



coherent behaviour below about 0.5 pc. As the findings of 
Goodman et al. are also based on the total linewidths this 
suggests that the change might rather reflect the transition 
from a regime where single separated clumps are identi- 
fied, to measurements of a superposition of substructures 
at smaller scales. 



3.2. Velocity probability distribution function 

Another quantity characterising the velocity structure 
both in observational data and in turbulence simulations 
is the probability distribution function (PDF) of veloci- 
ties. Although it contains no information on the spatial 
correlation in velocity space like the size-linewidth rela- 
tion or the A-variance, it shows complementary prop- 
erties, like the degree of intermittency in the turbulent 
structure (Falgaronc & Phillips 1990). The shape of the 
wings of the velocity PDF is thought to be diagnostic 
of intermittency, where the increasing degrees of inter- 
mittency produces a transition from Gaussian to expo- 
nential win gs. Two-dimensional Burg ers turbulence sim- 
ulations by Chappell fc Scalo (1999) , neglecting pressure 
forces, showed Gaussian velocity PDFs for models of de- 
caying turbulence and exponential wings for models driven 
by strong stellar winds. 

Due to the limited amount of information available 
from molecular lines, there is no direct way to deduce the 
velocity PDF from observations. One approach to deduc- 
ing the velocity PDFs is computation of the distribution of 
line centroid velocities (Kleiner & Dickmann 1985; Miesch 
& Bally 1994; Miesch et al. 1999). This method can also in- 
clude some information on spatial correlation as discussed 
in Sect. 3.3. However, the higher moments of the centroid 
PDF are very sensitive to the observational restrictions 
discussed in Sect. 2.2. 

Another method was introduced by Falgarone & 
Phillips (1990), who estimated velocity PDFs from 



high signal-to-noise observations of single line profiles. 
Investigating the statistical moments of profiles, Falgarone 
et al. (1994) found no simple Gaussian behaviour for 
many observations and provided a first comparison with 
three-dimensional (3D) hydrodynamic simulations. Most 
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of their PDFs could be represented by a superposition of 
two Gaussians where the wing component had about three 
times the width of the core component. Unfortunately, 
their method is only reliable for optically thin transitions 
at a very high signal to noise. We test both methods here, 
starting with the centroid velocity PDF. 

3.2.1. Centroid velocity PDFs 

In computing the centroid velocity PDF for a map one 
can either assign the same weight to each point in the 
map, or weight the different contributions by the inten- 
sities measured at that point. We find that the PDFs 
retain similar shape and the same wing behaviour with 
both methods, and therefore use intensity weighting in 
the following analysis, as it is less influenced by observa- 
tional noise. We have also used normal histograms here, 
instead of the more sophisticated Johnson PDF estima- 



10.000 r 



tor applied by Miesch et al. (1999) because the error bars 
present from the uncertainty about the noise treatment 
greatly exceed the influence of the numerical PDF esti- 
mator. 

Fig. U shows the centroid velocity PDFs for the three 
data sets. We find that the IRAM and the CfA data are 
characterised by an asymmetry of the velocity distribu- 
tion, indicating some kind of large-scale flow within the 
mapped region. Looking at the wings of the distributions, 
however, all three data sets are consistent with a Gaussian, 
which would appear as a parabola in the lin-log plots 
shown. Only at the scale of the CfA map is a definite 
conclusion not possible, due to the large error bars. 

Beyond this phenomenological approach, the shape of 
the PDFs can be quantified by their statistical moments. 
The most frequently used moments are 



(»c> 



K = 



dv c P(v c )v c 

) 

dv c P(v c )[v c - (v c )] 2 

> 

1 [°° 

-j / dv c P(v c )[v c - (v c )f 



(2) 
(3) 
(4) 



where (v c ) is the mean, a 1 the variance, and K the kurtosis 
of the distribution. The probability distribution function 
P is normalised to unity. The variance is a measure for 
the total turbulent mixing energy, while the kurtosis char- 
acterises the deviation from a Gaussian profile. It takes a 
value of three for a Gaussian distribution, and six for a dis- 
tribution with exponential wings. Values between about 
2.7 and 3.0 are still consistent with a Gaussian which is 
truncated due to a finite sample size. Shallow wings be- 
come obvious for kurtosis values above about 3.4. 

We restrict the analysis to these low moments because 
higher moments become increasingly uncertain due to the 
influence of observational noise, non-perfect spectral base- 
lines, and error-beam pickup (see e.g. Bensch et al. 2001b). 
Without extremely high signal to noise data it is impossi- 
ble to obtain reliable constraints on the spatial variation 
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Fig. 2. Probability density distribution of centroid veloci- 
ties for the three Polaris Flare data cubes. The error bars 
show the deviation introduced by different treatments of 
the observational noise, as discussed in the text. 



even for the next higher velocity moments. The situation is 
even worse for methods like t he spectral correlation func- 
tion ( Rosolowsky et al. 1999 ) that include all details of 
the line profiles. 

In Table [l] we give the standard deviation (square root 
of the variance) and kurtosis for the three PDFs. Only 
the KOSMA map shows an almost Gaussian distribution 
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Table 1. Parameters of the ccntroid velocity PDFs and the average line profiles in the three Polaris Flare maps 



PDF type 


telescope 


std. deviation 


kurtosis 


core width a 


wing width a 






[km/s] 




[km/s] 


[km/s] 


centroid PDFs 


IRAM 


0.17 ±0.01 


5.1 ±0.2 


0.130 ±0.005 


0.24 ± 0.05 




KOSMA 


0.53 ±0.03 


2.7 ±0.1 


0.57 ±0.05 


0.41 ±0.03 




CfA 


1.6 ±0.2 


5.0 ±0.3 


1.02 ±0.05 


2.2 ± 0.5 


line profiles 


IRAM 


0.88 


2.4 


0.98 ±0.04 


0.80 ±0.07 




KOSMA 


1.06 


2.7 


1.10 ±0.01 


1.03 ±0.05 




CfA 


2.3 


3.8 


1.8 ±0.1 


2.7 ±0.3 



a standard deviation of the Gaussian fit 



of velocity centroids with slightly steeper-than-Gaussian 
wings. The kurtosis for the other two maps is clearly larger 
than Gaussian. We have tested whether we can reproduce 
this by the superposition of two Gaussians as proposed 
by Falgarone et al. (1991). Within the error bars we al- 
ways obtained good fits. All wings can be reproduced by a 
Gaussian. The widths of the core and the wing component 
obtained from this fit are also given in Table @. In contrast 
to the width ratio of about three obtained by Falgarone 
et al. we find a ratio of about two. 

Only a few of the observed maps analysed by 



Miesch et al. (1999) showed approximately Gaussian cen- 
troid velocity PDFs, while the majority had PDFs with 
shallower wings that could be fitted with either ex- 
ponential laws or power laws. Fitting the exponent, 
Miesch et al. (1999) obtained two different values when 



treating cither the whole distribution or only the wings. 
This is similar to our fit of the observational data with two 
different Gaussians for the wing and the core component. 

In the Polaris Flare maps the wings of all distributions 
can be represented by Gaussians. However, the total dis- 
tributions sometimes deviate considerably from Gaussian 
behaviour resulting in a Gaussian kurtosis value at inter- 
mediate scales (the KOSMA map) compared to signifi- 
cantly larger kurtosis at smaller and larger scales. Besides 
the effect of observational errors, this might represent the 
influence of systematic velocity trends across the mapped 
region which appear mainly in the core component but not 
in the wing. We don't expect these trends in the isotropic 
simulations analysed in Sect. 



5.2 



3.2.2. Average line profiles 

In Fig. U we show the average line profiles for all three 
data sets. The corresponding moments and fit parameters 
of the distributions arc also given in Table |l|. The average 
line profile automatically contains the weighting of each 
velocity contribution by its intensity as discussed above 
for the centroid PDFs, so that we can compare both. 

The general shape is similar but the line profiles are 
much broader and the peak positions are not at the same 
velocity as in the PDFs. For a detailed comparison, we plot 
both the average line profile and the centroid velocity dis- 
tribution for the IRAM data on the same logarithmic scale 
in Fig. The ratio between average line profile width and 
centroid velocity PDF width measured from the Gaussian 









KOSMA 




CfA ; 







-2 
[km/s] 



Fig. 3. Average line profiles for the three Polaris Flare 
data cubes. The solid line shows the IRAM data, the dot- 
ted line the average profile in the KOSMA map multiplied 
by 2 and the dashed line the CfA data multiplied by 8. 



centroid velocity PDF 



average line prof 




Fig. 4. Comparison of the total velocity PDF given by the 
average line profile and the centroid velocity PDF for the 
IRAM data. 



fitted to the wings of the distributions is about 3.4 for the 
IRAM map, 2.2 for the KOSMA map, and 1.5 for the CfA 
map. Using the variance of the full distributions we get 
ratios of 5.2, 2.0, and 1.4 respectively for the three maps. 
The difference of the two ratios for the IRAM map cor- 
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responds to some large scale velocity flow on that scale 
producing the irregular PDF core seen in Fig. ^. 

The variation of this ratio with the size of the map is 
again naturally explained by the two different length scales 
involved: the line-of-sight integration and the size of the 
map. In the IRAM data, the small size of the map provides 
a relatively narrow centroid velocity PDF compared to the 
broad average line profile determined by the line-of-sight 
integration through the full depth of the cloud. In con- 
trast, the thickness of the cloud will certainly be smaller 
than the full extent of the CfA map. This is also indicated 
by the approaching slopes of the two different variances 
within beams of varying size in Fig. |] at large scales. To 
compare the ratios obtained here with turbulence simu- 
lations in model cubes we have to consider scales where 
the map size is about equal to the thickness of the cloud. 
From the intensity maps of the clouds and Fig. [j] we es- 
timate a thickness corresponding to about 2° in angular 
scale. The resulting typical value for the width ratio that 
we should reproduce in the turbulence simulations then 
falls between about 1.5 and 1.6. We will see that several 
but not all turbulence models show such values. 

We can resolve the long lasting dispute over whether 
to use the average line profiles or the centroid velocity 
PDF as a measure for the 3D velocity PDF. The answer 
is determined by the size scales involved in the observa- 
tions. Since the velocity centroids ignore the integration 
along the line of sight, they provide the correct distribu- 
tion when the map size is larger than or comparable to 
the thickness of the cloud, while for small maps, the line 
profiles provide the better average, because they include 
a larger sample from the longer line-of-sight integration. 

We have to mention two caveat s. First, optical depth 
effects can broaden lines. However, Bensch et al. (2001a) 



have found that the 13 CO and 12 CO Polaris data show 
the same spatial scaling laws, although the maps differ, 
suggesting that optical depth effects play only a minor 
role. Second, the two methods can only be equivalent for 
an isotropic medium, which is not guaranteed. 

3.3. Velocity difference PDFs 

From the velocity centroid maps one can also extract in- 
formation on the distribution of scales in the velocity field 
by considering PDFs of velocity differences between points 
separated by different lags (distances). This provides inde- 
pendent information on the structure of the velocity field. 
Investigation of the PDF of velocity differences as a func- 
tion of spatial separation has been pursued by Miesch & 
Scalo (1995), Lis et al. (1998), and Miesch et al. (1999). 
For a discussion of the details and the application to sev- 
eral molecular clouds we refer to Miesch et al. (1999). 

Here we don't study the full PDF of centroid velocity 
differences but the variation of the first statistical mo- 
ments of this PDF as a function of lag between the two 
points considered. Because of the symmetry of the velocity 
differences, all odd moments vanish. The first two non- 
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Fig. 5. The standard deviation of the velocity difference 
PDF as a function of lag. 

zero moments of the velocity difference distribution are 
the variance and the kurtosis: 



a 2 (L) 



K(L) = 



d 2 r / d 2 r' f(r)f(r')[v c (r)-v c (r')] 2 

map J \ r—r' \ — L 

d 2 r f d 2 r' f(r)f(r') 

map J \ r—r'\=L 

d 2 r f d 2 r'f(r)f(r')[v c (r)~v c (r')} 4 

map J \ r—r'\=L 

a\L) f d 2 r f d 2 r' f(r)f(r') 

J map J \r — r'\— L 



(5) 



(6) 



Integrations over the spatial vectors r and r' scan the 
whole map. The contribution of the velocity difference be- 
tween the points r and r' is weighted by weighting factors 
/. We compared equal weighting as used by Miesch et al. 
(1999), weighting by the geometric mean, and by the prod- 
uct of the two intensities, and find little difference. In the 
following we use the weighting by the geometric mean, as 
it is a linear intensity weighting for each term, as in the 
case of the PDFs discussed above. 

The quantity given by the variance of the two-point 
PDF a 2 (L) in Eq. (|) as a function of the lag between 
the points L is identical to the ordinary structure func- 
tion as used e.g. by [Miesch fc Bally (1994) , except for the 
normalisation of the structure function, so we can com- 
pare the results. Miesch et al. (1999) obtained for several 
clouds a power law behaviour for the variance of the cen- 
troid velocity differences cr 2 {L) cx L 7 , with 7 w 0.85 (0.33 
...1.05) except for the largest lags, where a 2 remained 
roughly constant. To enable a better comparison with the 
size-linewidth relation discussed above, we use here the 
standard deviation a instead of the variance . In Fig. I 
we show the resulting plot for the Polaris Flare data. We 
find an overall slope of 0.47, quite close to that found for 
the size-linewidth relation. The error bars are somewhat 
smaller but the two data sets at higher resolution show 
some decrease of the slope at the largest lag. This must 
be an artifact due to the finite map size, since it does not 
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Fig. 6. The kurtosis of the centroid velocity difference 
PDF as a function of the lag between the points consid- 
ered. 



continue at th e next larger scale. T he structure functions 
calculated by Miesch et al. (1999) show a much stronger 
flattening at large lags, going to constant values for all 
maps. This is probably due to the artificial removal of ve- 
locity structure at large lags introduced by their method 
to subtract large scale trends. 

The good agreement between the size-l inew idth re- 
lation for the centroid velocities in Sect. and the 
structure function discussed here seems inevitable when 
we consider that the variance in velocity differences on a 
certain scale is a kind of differential measure for the to- 
tal variance within a certain radius as measured with the 
scanning-beam size-linewidth relation. We thus expect a 
similar behaviour. 

The kurtosis K(L) is a measure for the correlation 
of the internal motions. Values exceeding three at dif- 
ferent lags indicate the stre ngth of the correlat ion in 
velocity space at those scales. Miesch et al. (1999) found 
that the velocity difference PDFs in the studied clouds 
change from kurtosis values between about 10 and 30 at 
small lags to nearly Gaussian behaviour at large lags. This 
is also typical for incompressible turbulence (She 1991). 
Lis et al. (1998) found strong non-Gaussian distributions 
at scales associated with filaments and approximately 
Gaussians at larger lags. 

In Fig. ^|the kurtosis of centroid velocity differences in 
the Polaris Flare data is plotted. In contrast to the other 
quantities, the different resolutions do not line up here to 
a single line but the kurtosis for each map drops inde- 
pendently to the Gaussian value at about the map size. 
This behaviour can be understood by considering which 
quantities at which scales determine the kurtosis. At the 
largest lag of any map the kurtosis measures mainly the 
shape of the centroid probability distribution of the whole 
map which is more or less cl ose to Gaussian for all data 
sets considered here (Sect. |3.2.2 ). Thus we can always 
expect a value around 3 when the scale for the kurto- 
sis determination approaches the map size. Contrary to 



the discussion provided by Miesch et al. (1999), this does 
not mean that there are no correlations at larger scales 
but that they cannot be addressed from points within the 
map. 

At all smaller lags the kurtosis is a measure for the 
correlated motions on that scale relative to the overall 
motions seen in the map, which is sca nned when comput- 
ing the kurtosis and variance. In Sect. 3J3 we will see that 
kurtosis values above three are produced only if the maps 
contain some motion on scales larger than the scale on 
which the kurtosis is measured. The steps in Fig. ^| are 
thus unavoidable when switching to another map since 
we always measure the correlated motions on a particular 
scale relative to the total motions in the map considered. A 
slightly sub-Gaussian behaviour at the largest scale might 
be produced by optical depth effects somewhat flattening 
the core of the distribution. 

3.4. The A-variance 

Stutzki et al. (1998) introduced the A-variance to measure 
the amount of structure present at different scales in multi- 
dimensional data sets. The A-variance at a given scale of 
an n-dimensional data set is computed by convolving the 
data with an n-dimensional spherical down-up-down func- 
tion of that scale, and measuring the remaining variance. 
The A-variance analysis computes the average variance on 
a certain scale similar to the structure function giving the 
variance of the velocity differences between two distinct 
points separated by a certain lag. For the A-variance, how- 
ever, the variance of the filtered map is computed, instead 
of the average variance of all point-to-point differences 
corresponding to a certain lag. Thus, the A-variance of 
a smooth map with a linear gradient vanishes, while the 
structure function discussed above detects the gradient. 
The advantage of the A-variance is its better sensitivity 
to specific spatial scales. It provides a good separation of 
systematic trends, structures on certain scales, and effects 
like noise. Furthermore, it allows the direct computation 
of the equivalent Fourier spectral index. A comprehensive 
discussion is given by Bensch et al. (2001a). A similar 
method was introduced recently by Brunt (1999) to char- 
acterise the 3D velocity structure of model cubes. He used 
a rectangular filter function composed of adjacent cubes of 
different size rather than the spherically symmetric filter 
function used to compute the A-variance. 

Bensch et al. (2001a) applied the A-variance analysis 
to the intensity maps of the Polaris Flare discussed above. 
We used the same method in paper I to analyse the den- 
sity structure in turbulence simulations, and compared the 
results to the observational data. Paper I also discussed 
the A-variance for the 3D velocity field of the simulated 
turbulence and compared it to the A-variance for the 3D 
density, but did not provide any direct comparison to the 
observations. We found that the A-variance of the veloc- 
ity behaves similarly to that of the density in showing 
the characteristic scale of the driving mechanism used in 
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Fig. 7. The square root of the A-variance of the velocity 
centroid maps of the Polaris Flare observations. In con- 
trast to paper I, we have chosen to take the square root, 
yielding the standard deviation of the filtered maps in- 
stead of the variance, in order to give a better comparison 
to the other linewidth related quantities. 

the turbulence models. However, the amount of structure 
observed at smaller scales differs between density and ve- 
locity. The turbulence creates many thin dense regions, 
leading to an exponent of the A-variance for the density of 
about 0.5, whereas it creates hardly any small-scale struc- 
tures in the velocity field, so that the A-variance for the 
velocity drops off much more steeply, with an exponent of 
about 2. 

Unfortunately this method of measuring the structure 
of the 3D velocity field in the simulations has no directly 
equivalent approach applicable to observations, since they 
only provide one-dimensional velocity information pro- 
jected onto the plane of the sky. There exists no simple 
relation between the three dimensional velocity structure 
and the behaviour of the projections. Therefore, we must 
instead apply the A-variance analysis to observable ve- 
locity parameters like the map of centroids which can be 
derived from both the simulations and the observations. 
They can be used to judge whether a simulation repro- 
duces observed properties. We do lose information by this 
procedure, of course. 

Fig. shows the square root of the A-variance for the 
Polaris Flare velocity centroid maps. The upturn at the 
smallest lags in each map is produced by the observational 
noise adding power at small scales which is only detected 
in the A-variance spectrum. Bensch et al. (2001a) showed 
how the influence of noise can be subtracted in the A- 
variance. If we apply this noise correction we find a turn- 
down at the smallest lags in each map. This lack of small 
scale variations in the different maps is due to the finite 
beam size of the observations blurring the structures. In 
the line intensity maps, the influence of the finite beam 
can be described by a convolution of t he structure with a 
Gaussian beam. Bensch et al. (2001a) provided a formal- 
ism to correct for both the observational noise and the 



beam smearing, if one assumes that the internal scaling of 
the structure is given by a power law. 

In the line centroid maps we cannot apply the same ap- 
proach because the beam convolution of the channel maps 
that is performed by a telescope does not correspond to 
the convolution of the line centroid map. Estimates for 
random velocity fields were obtained in simulations by 



Miesch & Bally (1994), but they do not provide a direct 
way to correct for the beam smearing in the observed data. 
To obtain an estimate for the uncertainty introduced by 
noise and beam smearing we have tested two extreme ap- 
proaches. If we apply the noise correction only and disre- 
gard the smallest lags in each map, which are obviously in- 
fluenced by the beam smearing, we obtain average slopes 
of 0.04 in the CfA map, 0.46 in the KOSMA map, and 
0.68 in the IRAM map. In this case the curves do not 
line up perfectly. The A-variance for one lag in a larger- 
scale map always falls somewhat below the A-variance for 
the same lag in the smaller-scale map, indicating that the 
slopes are probably overestimated. If, on the other hand, 
we assume that the beam provided a simple convolution 
of the centroid map, so that we can apply Eq. (15) from 



Bensch et al. (2001a), we obtain slopes of -0.05, 0.36, and 



0.32, respectively and the opposite misalignment. The true 
slopes should fall in between these two extreme cases. 

Over most scales covered by the CfA map, the slope 
flattens to zero, in contrast to the behaviour in the inten- 
sity maps. This virtual lack of large-scale variations does 
not reflect the real structure but appears to be due to 
the missing weighting in the A-variance analysis. Here, 
all points are counted equally, even those with intensities 
well below the noise limit. Thus, the A-variance analysis 
necessarily fails in cases where the maps are only sparsely 
filled by emission. Less than one third of the CfA map 
shows emission above the noise limit so that we cannot 
expect any significant results at scales above half a de- 
gree. 

For the maps at smaller scales without large "empty" 
regions, we can expect significant A-variance spectra but 
we must admit that without a good theory for the in- 
fluence of the beam convolution on line centroids we can 
hardly separate the effects of the beam smearing from the 
velocity scaling in maps which are not large compared to 
the beam size. Thus large error bars remain for the mea- 
sured velocity scaling in the Polaris Flare preventing any 
definite conclusion on the scale dependence of the slope in 
the velocity structure. 

In contrast the size-linewidth relation and the struc- 
ture function do not show clear signatures of noise and 
beam smearing in Figs, [l] and g. To test their sensitivity 
we have applied them to data smoothed with a Gaussian 
filter tuned to ensure that the corresponding A-variance 
shows the same behaviour a s the A-variance whe re the 
noise was removed following Bensch et al. (2001a) . Only 



small changes relative to the original data are seen. The 
size-linewidth relation and the structure function of the 
velocity centroids exhibit a weak steepening of the slope 
by about 0.04, while the velocity PDFs and the kurtosis 
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show no clear differences. The comparison of the results 
from the original maps and the smoothed data shows that 
these other methods are not particularly sensitive to un- 
correlated noise or beam smearing at small scales. 



3.5. Comparison of the methods 

We can classify all of the methods we have described 
in terms of the velocity information setup in a two- 
dimensional map, the filtering function used, and the 
weighting of the data in the map. 

Most analyses were restricted to the velocity centroids, 
which are effectively the first moment of the local veloc- 
ity profile. The size-linewidth relation adds the local vari- 
ance, i.e. the second velocity moment, and the study of 
the PDFs also uses the kurtosis. However, with sufficiently 
high signal-to-noise, higher moments may provide valuable 
additional information. 

The variance is strongly dominated by the depth of the 
observed cloud, so that it contains information lost when 
considering the velocity centroids only. We have seen that, 
for maps where the line profiles sample the cloud deeply in 
comparison to the map size, the integrated line profile is a 
better measure for the true velocity distribution function 
than the PDF of the velocity centroids. 

We have applied three different kinds of filters: 
the scanning-beam size-linewidth relation effectively con- 
volves the map with a positive Gaussian filter; the A- 
variance analysis uses a spherically symmetric up-down 
filter; and the structure function uses a filter consisting of 
a positive and a negative spike separated by a certain dis- 
tance. In the latter case, spherical symmetry is provided 
by the superposition of the resulting variance values for 
different directions of the filter axis. The structure func- 
tion is sensitive to large-scale gradients and can detect 
certain geometric structures, but because of the strong 
localisation of the filter in the spatial domain, it is unfor- 
tunately sensitive to a broad spectrum of spatial frequen- 
cies in the Fourier domain. In the statistical analysis of 
velocity fluctuations, it is therefore at a disadvantage in 
the detection of characteristic scales and frequencies com- 
pared to the A- variance analysis. Similar conclusions were 



obtained by Houlahan & Scalo (1990). This may explain 



why we detect the scale of noise only by means of the A- 
variance, whereas its influence is hidden in the full spectra 
of the structure function and the size-linewidth relation. 

The A-variance analysis, on the other hand, does not 
yet take into account different weights for the informa- 
tion in different regions of an observed map, so that it 
fails for maps with large regions dominated by noise. The 
weighting of the velocity centroid information by the in- 
tensity, as is done automatically in the size-linewidth re- 
lation, reduces the uncertainty due to observational noise 
when computing the centroid probability distribution or 
structure function. 



3.6. Other approaches 

With sufficiently high signal-to-noise ratio, it is possible 
to extend the methods discussed here. Overviews of the 
different existing methods have been provided recently 
by Vazquez-Semadeni (2000) and Ossenkopf et al. (2000). 
First, one can apply the basically two-dimensional meth- 
ods to higher moments of the line profiles, providing new 
information especially on the intermittency in velocity 
space. Alternatively, the velocity channel maps can be 
analysed as demonstrated with the A-variance analysis 
by Ossenkopf et al. (1998). Another method is to com- 
pare full spectra, using the spectral correlation function 
(Rosolowsky et al. 1999) and extending this method to 
consider all spatial variations. 

Tauber (1996) discussed the smoothness of line profiles 
as a measure for the size and number of coherent units con- 
tributing to the profiles. Applying a rough approximation 
to this analysis, Falgarone et al. (1998) conclude that the 
size of cells in the Polaris Flare observations must be as 
low as 200 AU. 

When looking for characteristic global features in the 
density-velocity structure, the principal component analy- 
sis introduced by Heyer & Schloerb (1997) is probably the 
most significant tool. It identifies the main components in 
the position-velocity space in terms of eigenvectors and 
eigenimages. Although the principal component analysis 
represents a reliable method to find the dominant struc- 
tures even in complicated images, the significance of the 
higher-order moments still has to be determined. 



4. Turbulence models 

4.1. Simulations 



We use simulations of uniform decaying or driven tur- 
bulence with and without magnetic fields described by 
Mac Low et al. (1998) in the decaying case and by Mac 
Low (1999) in the driven case. These simulations were 
performed wi th the astrophysical MHD code ZEUS-3D0 
( Clarke 1994 ). This is a 3D version of the code described 
by Ston e fc Norman (1992a, b) using second-order ad- 
vection ( van Leer 1977 ), th at evolves magnetic fie lds us- 
ing constrained transport ( Evans fc Hawley 1988|) , mod- 
ified by upwinding along shear Alfven characteristics 
( Hawley fc Stone 1995 ). The code uses a von Neumann 
artificial viscosity to spread shocks out to thicknesses of 
three or four zones in order to prevent numerical instabil- 
ity, but contains no other explicit dissipation or resistivity. 
Structures with sizes close to the grid resolution are sub- 
ject to the usual numerical dissipation, however. In Paper 
I we discussed the effects of limited numerical resolution, 
which leads to numerical viscosity, and noted that resolu- 
tion studies could be used to determine which properties 
were well resolved. 



2 Available from the Laboratory for Computational 
Astrophysics of the National Center for Supercomputin g 



Applications, http://zeus.ncsa.uiuc.edu/lca_home_page.htm] 
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The simulations used here were performed on a 3D, 
uniform, Cartesian grid with side L = 2 and periodic 
boundary conditions in every direction, using an isother- 
mal equation of state. To deal with velocities compara- 
ble to those in the observations we have assumed here a 
cloud temperature of 10 K corresponding to a translation 
of the dimensionless sound speed in the simulations to a 
physical sound speed of 0.2 km/s in the data analysis. 
The initial density and, in relevant cases, magnetic field 
are both initialised uniformly on the grid, with the initial 
density po = 1 and the initial field parallel to the z-axis. 
The turbulent flow is initialised with velocity perturba- 
tions drawn from a Gaussian random field determined by 
its power distribution in Fourier space, as described by 
Mac Low et al. (1998). For decaying models we use a flat 
spectrum with power in the range fed = 1 — 8, where the 
dimensionless wavenumber fed = L/X c i counts the number 
of driving wavelengths in the box. For our driven mod- 
els, we use a spectrum consisting of a narrow band of wave 
numbers around some value k^, and driven with a fixed 
pattern at constant kinetic energy input rate, as described 
by Mac Low (1999). 

We have tested the influence of numerical viscosity by 
running the simulations with the same physical parame- 
ters on grids of 64 3 , 128 3 or 256 3 zones. Higher resolution 
grids have numerical viscosity acting at smaller scales, so 
changing the resolution shows the effects of numerical vis- 
cosity on our results. The influence of numerical resolution 
on the simulation results is discussed below separately for 
each for the statistical measures. 

4.2. Simulated Observations 

To compare our simulations with observations we must 
synthesise observational maps from the simulated density 
and velocity fields. We assume that the cubes are opti- 
cally thin for this first study so that direct integration 
along lines of sight through the cube neglecting optical 
depth effects yields line profiles. This appears to be a rea- 
sonable assumption for comparison with the low column 
density clouds in the Polaris flare observed in 13 CO, but 
is a worse assumption for higher column density clouds or 
more optically thick species. 

We also neglect the periodic nature of the simula- 
tions, effectively observing the simulation cubes as isolated 
structures in a vacuum. This second assumption must be 
taken into account in analyses affected by the path-length 
through the cloud, such as comparisons of the velocity 
PDF measured from the average line profile vs. the cen- 
troid velocity distribution. 

4.3. Statistical fluctuations 

To get a feeling for the significance of the structural prop- 
erties indicated by the different measures relative to the 
statistical variations in the turbulence, which may build 
up from isotropic initial perturbations during the turbu- 
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Fig. 8. A-variances (upper panel) and structure functions 
(lower panel) in the ccntroid velocity maps obtained for 
a decaying hydrodynamic model (model D from Mac Low 
et al. 1998) in the three possible directions of projection. 
The dotted line indicates a slope of 0.5 for comparison. 
The variation between the different directions is similar 
for the other models. In the A-variance plots we always 
give the standard deviation (the square root of the A- 
variance) for a better comparison to the other velocity 
variations. 

lence evolution, we can compare different directions within 
the same model cube. 

Fig. H shows the A-variances and structure functions 
for the three centroid velocity maps along the three axes 
in a decaying hydrodynamic model. The A-variances show 
variations of 40% in value at the largest scales, where the 
exact modal structure is still different in the three direc- 
tions. The structure function shows smaller variations cor- 
responding to its lower sensitivity to changes at particular 
spatial scales. The variation of the slope in the different 
directions falls below 0.1. The equivalent plot for the size- 
lincwidth relation is similar to that for the structure func- 
tion. The variations in the kurtosis, on the other hand, 
are more like those seen in the A-variance. Here, we some- 
times find small dips and rises distorting the monotonic 
decay to values slightly below 3 at the largest scales. In the 
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Fig. 9. Size-linewidth relation for the velocity centroids 
of a hydrodynamic simulation driven at wavenumbers of 
fed = 2 (solid), kd = 4 (dashed), and fed = 8 (dot- 
dashed). These are models HE2, HE4, and HE8 from 
Mac Low (1999). The dotted line shows a power-law with 
slope 0.5. 



Fig. 10. Size-linewidth relation for the velocity centroids 
in a hydrodynamic simulation of decaying turbulence at 
time t — 0.75t C ross (Across is the initial crossing time of the 
region), computed at resolutions of 64 3 (dot-dashed), 128 3 
(dashed), and 256 3 (solid). These are models B, C, and D 
from |Mac Low et al. (1998)| . 



total and the centroid velocity PDFs we find substantial 
variations in the core and the central position of the dis- 
tributions corresponding to the different largest velocity 
modes but no changes in the wing behaviour. 

5. Statistical description of simulations and 
comparison to observations 

5.1. Size-linewidth relation 

We begin by considering the results of applying the size- 
linewidth analysis described in Sect. 3.1 to the models. In 
Fig. H we show the size-linewidth relation for the velocity 
centroids in three models of hydrodynamic driven turbu- 
lence that differ only in the scale that they are driven. The 
driving wavelengths are 1/2, 1/4, and 1/8, respectively. 

We find power-law behaviour through most of the 
regime only for the model driven at the largest available 
scales. Models driven with smaller characteristic scales 
show a flattening of the relation at lags above the driving 
scale. A slight flattening at the largest lags is also visible 
in the observational data. This appears to be an indica- 
tion of a turbulence driving scale close to the size of the 
molecular cloud. 

A drop off in velocity dispersion is seen at small lags 
in all of the models. This can be explained straightfor- 
wardly as an effect of numerical viscosity. In Fig. [l^ we 
show a comparison of three models of decaying turbulence 
that are statistically identical, but were computed at res- 
olutions of 64 3 , 128 3 , and 256 3 zones. Increasing resolu- 
tion results in decreasing numerical viscosity, so Fig. [To] 
demonstrates explicitly the effect of changing the numeri- 
cal viscosity. As can be seen, the slope does not change at 
large lags, and the higher resolution models agree within 
a few percent on the magnitude of the velocity dispersion. 



At small lags, on the other hand, the velocity dispersion 
falls off. This occurs in all models at roughly the same 
number of grid zones, and hence at larger physical scale 
in the lower resolution models. This part of the spectrum 
thus reflects the effect of numerical diffusion eliminating 
some of the small-scale structure. The same drop-off oc- 
curs in each model (see Fig. ||), so the behaviour at the 
smallest scales should be viewed with caution. We note, 
however, that a physical diffusion such as ambipolar dif- 
fusion (Zweibel & Josafatsson 1983, Klessen et al. 2000) 
is expected to produce similar behaviour at the diffusion 
scale. 

Magnetic fields do not appear to modify the size- 
linewidth relation, although they can make order unity 
differences in the magnitude of the velocity dispersion and 
produce significant anisotropy, as described below in Sect. 

n 

We can conclude that the observed power-law be- 
haviour of the size-linewidth relation is reasonably ex- 
plained by either hydrodynamic or magnetised turbulence 
driven at scales comparable to the largest observed scales. 

5.2. Velocity probability distribution function 
5.2.1. Centroid velocity PDFs 

We begin our study of velocity PDFs in the models by 
examining the effect of numerical resolution on the cen- 
troid velocity PDFs. In Fig. |ll| the 3D PDFs are com- 
pared to the centroid PDFs for a model of decaying MHD 
turbulence at different resolutions. All PDFs are well rep- 
resented by Gaussians. The width of the PDF drops at 
smaller resolution. This can be attributed to the stronger 
influence of the numerical viscosity in the lower resolution 
cubes damping the turbulence faster. This agrees with the 
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Fig. 11. (a) 3D and (b) centroid velocity PDFs in a model 
of strongly magnetised decaying turbulence with initial 
Mach number M — 5 and Alfven number A = 1 after 
1.5t cro ss at a resolution of 256 3 (solid), 128 3 (dashed), and 
64 3 (dashed-dotted). (Models Q, P, N from Mac Low et 
al. 1998) 



measurements of resolution effects on kinetic energy de- 
scribed in Mac Low et al. (1998), in which the magnitude 
of the kinetic energy increased with resolution, although 
the decay rate was constant. An additional effect is seen 
in the centroid PDFs. At low resolutions the sampling of 
the wings of the Gaussian is insufficient, so that the dis- 
tribution appears too narrow. The kurtosis values of these 
PDFs are 2.6, 2.8, and 2.9 with growing resolution. Hence, 
sub-Gaussian kurtosis values can at least partially be ex- 
plained by small map sizes. 

In Table | we give the PDF moments for most of 
the models discussed, covering a wide range of different 
physical properties. The first several columns describe the 
model input parameters, and the remaining columns con- 
tain the parameters of the PDFs obtained. 

We first consider the effect of varying the driving 
wavenumber, holding the energy input constant. Fig. [12| 
and the corresponding values in Table || show that driv- 
ing at smaller wavenumbers (longer wavelengths) produces 



Fig. 12. Centroid velocity PDFs from three models driven 
with the same strength at wavenumbers of fed — 2 (solid) , 
fed = 4 (dashed), and fed = 8 (dash-dotted) (models HC2, 
HC4, and HC8 from Mac Low 1999) 



broader PDFs, because such models have lower dissipa- 
tion rates (Mac Low 1999), and thus higher rms velocities. 
More interestingly, we find that driving at the largest scale 
in the model fed = 2 produces a centroid velocity PDF 
with apparently non-Gaussian shape, as reflected in the 
kurtosis value of 2.4 given in Table A similar result was 
found by Klessen (2000), who argued that it is most likely 
due to cosmic variance. That is, an insufficient number 
of modes arc sampled at these long driving wavelengths 
to fully describe a Gaussian field, so the PDF appears to 
have a distorted shape. Depending on the random num- 
bers used to initialise the largest modes, both Gaussian 
and non-Gaussian kurtosis values are then possible. 

To demonstrate the effect of magnetic fields, Fig. [l| 
shows the centroid velocity PDFs from a model of de- 
caying magnetised turbulence. It clearly indicates an 
anisotropic decay, with velocity components perpendicular 
to the magnetic field decaying substantially more quickly 
than velocities parallel to the field. In both cases, though, 
the PDFs remain Gaussian even at late times. These con- 
clusions are quantitatively supported by Table |[ 

We also examined magnetised driven turbulence. The 
fields do slightly shift the peak, but, as shown in Table |[ 
they still produce only marginally non-Gaussian centroid 
velocity PDFs in cases where the hydrodynamic model is 
Gaussian, contrary to the speculation of Klessen (2000) 
that magnetic fields might be an important alternative 
cause of non-Gaussian PDFs. The PDFs observed parallel 
to the field are roughly 20% wider than the perpendicular 
observations as seen in Table § 

Klessen (2000) showed that driving from large scales 
can produce non-Gaussian PDFs, but worried that every 
additional piece of physics also appeared likely to produce 
non-Gaussian PDFs, allowing no conclusions to be drawn 
from their occurrence. We have, however, demonstrated 
that neither magnetic fields nor the vorticity introduced 
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Table 2. Parameters of the ccntroid velocity PDFs and the total velocity PDFs for the model cubes 
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a Driven models (H & M series) from Mac Low (1999)| ; decaying models (single letters) from 
b Mechanical driving luminosity in arbitrary units (see Sect. 2.3 of Mac Low (1999) for all unit conversions) 
c Driving wavenumber. 
d Number of zones in each dimension. 

6 rms Mach number: initial value for decaying models, equilibrium value for driven models. 
f Ratio of Alfven velocity to sound speed. 

g Time at which values are measured, in sound-crossing times. 
h Variance of distribution for line centroid velocities and full cube. 
J Kurtosis of distribution for line centroid velocities and full cube. 
k Unpublished model. 



by shock interactions in driven or decaying turbulence pro- 
duce strongly non-Gaussian PDFs. 

Therefore, the non -Gaussian PDFs observed in the 
Polaris Flare (Sect. |3.2| ) must have a different explanation. 
Another candidate for producing non-Gaussian PDFs is 
self-gravity, but the lack of star-forming activity in the 



Polaris Flare suggests that self-gravitation does not play a 
dominant role there. This suggests that the non-Gaussian 
PDFs observed there are indeed due to the cosmic vari- 
ance introduced by driving from the largest scales of the 
region. The driving scale may actually be even larger than 
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Fig. 13. Centroid velocity PDFs from a model of decay- 
ing, magnetised turbulence (model Q from Mac Low et al. 
1998) at times in units of the initial crossing time t CTOSS of 
0.5 (solid), 0.75 (dashed), 1.5 (dash-dot), and 2.5 (dotted) 
observed (a) perpendicular to the field and (b) parallel to 
the field. 

identified here because we have not used any information 
from the atomic gas at larger scales. 

5.2.2. Average line profiles 

Using the turbulence simulations we can now revisit the 
question from Sect. 3.2.2 of how best to measure the ac- 
tual 3D velocity PDF from observations. In Fig. |l4| we 
compare the PDF of the whole velocity distribution that 
would be measured as the average line profile in an op- 
tically thin LTE medium to the PDF of the centroid ve- 
locities for a hydrodynamic model driven at fed = 2. This 
plot may be compared with Fig. ||, showing the centroid 
velocity distribution and the average line profile measured 
for the IRAM map in the Polaris Flare. 

In both cases we find similar distorted Gaussian distri- 
butions, with the centroid velocity distribution narrower 
than the full velocity distribution. Table ^ contains the ra- 
tio between the widths of the two distributions for all mod- 
els. The model ratios show substantial variation around a 
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Fig. 14. Comparison of the total velocity PDF and the 
centroid velocity PDF for a hydrodynamic model driven 
with one tenth of the energy of the example above (model 
HC2 from Mac Low 1999). 



typical value of about two. The strongest systematic varia- 
tion appears to be with initial or driven rms Mach number 
Mo, with higher Mq giving lower ratios, down to as low 
as 1.1 for the most strongly driven model HE2. Observing 
parallel to the field lines leads to somewhat lower ratios 
than perpendicular, presumably due to the higher veloc- 
ity variance seen along parallel lines of sight. However, 
the ratio also decreases during the decay of turbulence for 
unclear reasons. 

The ratio between depth and width is fixed to unity 
in the model. From the combination of the three Polaris 
Flare maps we estimate in Sect. 



3.2.2 a typical width ratio 



between 1.5 and 1.6 for a depth comparable to the lateral 
extension. In Fig. |l4|, showing a very high Mo simulation, 
the ratio is 1.5, but for the majority of simulations we find 
ratios of over 1.8, suggesting that the Polaris Flare obser- 
vations are of a region containing hypersonic turbulence 
with Mo of order 10 or higher. 

For most model PDFs we have tested three possi- 
ble fits to the distribution: Gaussians fitting either the 
whole distribution or only the wings, and an exponential 
fit. Although exponential wings cannot be ruled out com- 
pletely, the Gaussian fits are clearly better. We have found 
this to be true for all models of driven and decaying tur- 
bulence discussed here, as can be seen by visual inspection 
of Figs, [ll] to [l4|. This result agrees with the Gaussian ve- 
locity PDFs found by Chappell & Scalo (1999) for models 
of decaying Burgers turbulence. 

5.3. Velocity difference PDFs 

5.3.1. Second moments: structure function 

In Sect. 3J3 we have shown that the second moment of 
the centroid velocity difference PDFs, or equivalently the 
structure function, is a differential measure for the size- 
lincwidth relation, so that the same power law was found 
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Fig. 15. Standard deviation of the centroid velocity differ- 
ence PDF as a function of lag, which is equivalent to the 
structure function, normalised to the square root of the 
lag for 3D velocity cubes of a hydrodynamic simulation 
driven at wavenumbers of fed = 2 (solid), fed = 4 (dashed), 
and fed = 8 (dot-dashed). These are models HE2, HE4, 



and HE8 from Mac Low (1999) 



for the two functions when applied to the observations. In 
Fig. [l5| we plot the structure function of the three models 
driven at different scales whose size-linewidth relation was 
shown in Fig. 

Comparison of the two figures demonstrates that the 
structure function provides a better estimate of the driv- 
ing scale as its curvature is restricted to a narrower range 
around the driving scale than the curvature of the size- 
linewidth relation. On the other hand, only the latter in- 
dicates the dissipation scale by curvature at small lags. 
The different sensitivity of the two functions is caused by 
the different shape of their effective filter functions (Sect. 
3.5). We expect that with sufficient dynamic range there 
would be a set of lags where both functions would show 
power-law behaviour, but our 128 3 models have such lim- 
ited dynamic range that there is effectively no scale where 
this is true. The combined set of observations, on the other 
hand, does have enough dynamic range to be dominated 
by the power-law slope, so that both functions agree in 
most parts of the spectrum. 

5.3.2. Fourth moments: kurtosis 

The kurtosis of the centroid velocity difference distribu- 
tion can measure correlations in the motion at certain 
scales. Values around three indicate Gaussian distribu- 
tions, implying motions uncorrelated relative to the overall 
velocity field of the map considered. Miesch et al. (1999) 
found a decay of the kurtosis from small to large scales 
roughly proportional to the square root of the scale down 
to a valu e of three at the largest scales in their maps. In 
Sect. 3.3 we found somewhat shallower slopes, and a flat 



Fig. 16. Kurtosis of the centroid velocity difference dis- 
tribution for a hydrodynamic simulation driven at scales 
of kd — 2 (solid), fed = 4 (dashed), and fed = 8 
(dash-dot). Thes e are models HE2, HE4, and HE8 from 
[Mac Low (1999) . 



Gaussian section at scales above the total cloud size in the 
CfA Polaris Flare map. 

In Fig. |l6| we test how different driving scales in a 
driven hydrodynamic model influence the resulting kurto- 
sis plots. We clearly see that the models driven at larger 
scales reach the Gaussian value of kurtosis at larger lags, 
about a factor of two below the peak size of the driving 
structure indicated by the A-variance. This suggests that 
motion at lags above the scale of the driving process re- 
mains uncorrelated. At smaller lags, our uniformly driven 
models show kurtosis following As -1 / 3 , similar to the be- 
haviour observed for the Polaris Flare maps in Sect. |3.3| . 
None of our models reached the kurtosis values around 50 
shown by the observations of Miesch et al. (1999), sug- 
gesting that additional physics, especially in the driving 
function, may be responsible for these high values. 

Varying the resolution of our simulations does not ap- 
pear to markedly change the peak value of kurtosis, al- 
though the scale at which that value is reached is always 
the smallest scale of the simulation, suggesting that the 
correlated motions introduced by dissipation are a major 
influence in producing non-Gaussian velocity-difference 
PDFs. Magnetic fields do change the shape of the ve- 
locity difference PDF, though not drastically. Increasing 
field strength increases the peak value of kurtosis by about 
30%. 

5.4. A-variance 

In paper I we used the A-variance to investigate the spatial 
scaling of turbulent density structure. Here, we use this 
method to investigate the velocity structure. In contrast 
to our analysis of the observations, we have access to the 
full 3D structure in the simulations. 
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Fig. 17. A-variance for the velocity centroids for 3D ve- 
locity cubes of a hydrodynamic simulation of decaying tur- 
bulence at time t = 0.75t croB8 , computed at resolutions of 
64 3 (dot-dashed), 128 3 (dash ed), and 256 3 (solid). T hese 
are models B, C, and D from Mac Low et ah (1998). 



Fig. 18. A-variances computed for 3D velocity cubes of a 
hydrodynamic simulation driven with wavenumber fed = 2 
(solid), fed = 4 (dashed), and fed = 8 (dot-dashed). These 



are models HE2, HE4, and HE8 from Mac Low (1999) 



We start by studying how changes in the resolution, 
and so in the dissipation length scale, appear in the A- 
variance spectrum. In Figure ^ we compare models of 
decaying turbulence that are statistically identical, but 
were computed with resolutions of 64 3 (dash-dotted), 128 3 
(dashed), and 256 3 (solid). The variations at large scales 
(As > 0.1) are most likely due to statistical fluctuations 
rather than the changes in resolution. The result of sys- 
tematically increasing the resolution can be seen at small 
scales, as the spectrum reaches smaller scales at higher res- 
olution. The shape of the spectrum at the smallest scales 
in each model is very similar, with only the scale chang- 
ing. This shows the range over which numerical diffusion 
is acting. Above that scale, the models agree fairly well in 
a region that can be considered to be the inertial range 
of the turbulence. We see that, equivalent to the density 
structure, numerical diffusion causes a steepening of the 
spectrum, and that this reaches to a scale of roughly ten 
zones in every model. 

Now we can consider the effects of different driving 
wavelengths. Fig. [l8| shows the results of the analysis for 
three hydrodynamic models driven at different scales (the 
A-variance for the density distribution in these models 
is shown in Fig. 7b of paper I). We plot the square root 
of the A-variance, as that is the quantity directly related 
to the linewidth. We see pronounced peaks close to the 
driving scale at 0.5-0.6 Ad- They are even somewhat more 
pronounced than in the density structure of these models 
(paper I), because the driving process itself is implemented 
in velocity space (Mac Low 1999). Below the driving scale 
we find power-law behaviour with a slope of 0.57 for Av 
corresponding to a slope of 1.14 for the A-variance, down 
to the scale of the numerical viscosity at about ten pixels 
(As « 0.04 in the figure). 



How does this slope compare to the expected power 
spectrum? As discussed in paper I, there is a theoretical 
relation between the slope of the A-variance a and the in- 
dex £ of an n-dimensional power spectrum P(k) oc |fc| -< = 
given by a = ( — n. Please note the difference to the 
index (i nt of the often used one-dimensional power spec- 
trum P(k) that is obtained by angular integration of P(k), 
C = Cint +n—l. Hence, the slope is translated into values 
of 4.14 and 2.14 for £ and (i n t, respectively. Incompressible 
Kolmogorov-type turbulence is characterised by values of 
11/3 and 5/3, respectively, while shock-dominated veloc- 
ity fields should show the Fourier transform of a step func- 
tion, i.e. values of 4.0 and 2.0, respectively. The results 
from the simulations are quite close to the behaviour of 
shock-dominated gas, but the small deviation might ask 
for further investigation. 

For the density structure, we found A-variance slopes 
between 0.45 and 0.75 in projection, corresponding to C 
(not £i nt ) values between 2.45 and 2.75, considerably flat- 
ter than the slope in the velocity structure found here. 
We speculate that this is due to the compressibility of the 
isothermal gas modelled here. The gas piles up in thin 
sheets and filaments, so that the density structure has be- 
haviour somewhere in between that of a S function and a 
step function. The Fourier transform of a one-dimensional 
5 function in the density cube is characterised by £ = 2. 
The velocity on the other hand, remains uniform for longer 
distances behind shocks, so that its steeper spectral slope 
approaches more closely the value expected for a box full 
of pure step functions. 

We can compare the 3D A-variance spectrum to that 
from simulated velocity centroid maps. In Fig. |l^ we show 
the A-variance plots for the centroid velocity maps from 
the models shown in Fig. [H| Here, identification of the 
driving scale is no longer so easy. Nevertheless, the gen- 
eral functional behaviour is similar. The broad maximum 
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Fig. 19. A-variances of the centroid velocity maps of the 
three models from Fig. |l8|, with driving wavenumber fed = 
2 (solid), fed = 4 (dashed), and k^ = 8 (dot-dashed). 



Fig. 20. A-variance of the 3D velocity structure in a de- 
caying MHD turbulence model at times of 0.5, 0.75, 1.5, 
and 2.5 initial crossing times. 



is shifted depending on the driving wavelength, and there 
is a power law range with the slope close to 0.5. Thus, 
the A-variance of the velocity centroids does seem to re- 
flect the true cloud velocity structure, if not as clearly 
as the 3D spectrum. In the centroid velocities we mea- 
sure the same drift index a as in the 3D data cube, in 
contrast to the density structure where the integrating 
projection to two dimensions preserves the power spectral 
ex ponent Q, and ther efore increases the drift index a by 
1 ( Stutzki ct al. 1998 ). This conclusion also holds for the 
other models so that we restrict ourselves here to the 3D 
A-variance, noting that the observational centroid maps 
do not reveal the structural properties as clearly as dis- 
cussed here. 

Comparing Fig. |l^ with Figs. || and 15 or Fig. [ll] with 
Fig. [h], respectively, we see that the A-variance of the ve- 
locity centroids resolves both the driving scale and the dis- 
sipation scale most clearly. Although the structure func- 
tion also shows the driving scale by a change of the slope 
from 0.5 to 0, the bending from 0.5 to about -0.5 in the 
A-variance (Fig. |l9|) can be detected more easily in noisy 
data. The steepening at the dissipation scale which is also 
visible in the size-linewidth relation by a change of the 
slope from 0.5 to about 0.8 is also greatest for the A- 
variance where a slope of about 1.2 is measured at the 
smallest scales. 

In Fig. |(] we show the A-variances for the 3D velocity 
structure in a sequence of four time steps in the evolution 
of strongly magnetised (ua/c s = 10), decaying turbulence. 
The magnitude of the variance drops over time, starting 
from small scales, resulting in an increase in the effective 
peak of the driving function and a slight increase in the ef- 
fective slope of the spectrum over time. This effect was also 
seen in the density structure in paper I, and there shown to 



be proportional to t 1 / 2 as predicted by Mac Low (1999) 



To demonstrate the effect of magnetic fields we show 
in Fig. pi] the 3D A-variance plots observed parallel and 
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Fig. 21. A-variance of the 3D velocity structure observed 
parallel (thick lines) and perpendicular (thin lines) to the 
magnetic field for a models driven with k^ = 4 and no field 
(solid), weak field (va/c s — 1, dashed), and moderately 
strong field (va/c s = 5, dash-dot). 



perpendicular to the initial magnetic field in two MHD 
simulations with different field strengths, compared to a 
corresponding hydrodynamic simulation. We find a sub- 
stantial anisotropy in the case of a strong magnetic field 
(as in Sect. 



5.2.l| ), in contrast to the hydrodynamic sim- 



ulation and weak field cases. This is visible in the figure 
which shows a velocity variance for the strong-field case 
(dash-dot curves) substantially greater when viewed par- 
allel to the initial field direction (thick line) than perpen- 
dicular to it (thin line). 

The theory of incom pr essible MHD turbulence 



Sridhar & Goldrcich 1994, Goldrcich & Sridhar 1995 



|Goldrcich fc Sridhar 1997 ) predicts that the perpendicu- 
lar cascade will be more efficient than the parallel cascade, 
draining energy more quickly from perpendicular motions 
and producing anisotropic wavevectors aligned parallel 
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to the field. Our results suggests that this continues to 
be valid in the compressible regime. Further support for 
this conclusion may be drawn from the elongation of the 
density structures also seen in the simulations (see Fig. 5 
in Mac Low 1999). 

In the weak field case (dashed line) the anisotropy is 
far smaller, and indeed reversed, with perpendicular mo- 
tions having slightly greater power than parallel motions. 
However, in this case, unlike the strong field case, the field 
is weak enough that it is carried by the flow into a tangled 
geometry, so that the initial field direction no longer de- 
scribes the geometry of the field well. The small anisotropy 
observed is only a stochastic effect of the particular field 
configuration at that time. In fact, the plotted model is 
the only example from our simulations where we find this 
reversal even at the low field strengths. 

Aside from the anisotropy, these driven models show 
little overall difference in the velocity structure between a 
hydrodynamic model and the corresponding MHD models. 
A small increase in the amount of short-wavelength struc- 
ture is seen with increasing magnetic field, as in the den- 
sity structure (paper I), although there is no substantial 
shift of the effective driving scale. The small magnitude of 
this effect supports the claim by Heitsch et al. (2000) that 
the transfer of energy to smaller scales is insufficient to 
support the smaller scales against gravitational collapse 
if they were able to collapse in the absence of magnetic 
fields. 



6. Physical parameters 

6.1. Observations 

The velocity structure in the Polaris Flare maps shows 
basically a power-law behaviour below a few parsecs. At 
the smallest scales we find some indication for a steepen- 
ing of the A-variance spectrum of the velocity centroids, 
but the noise level in the data prevents a definite conclu- 
sion here. A similar scaling behaviour was demonstrated 



by Bensch et al. (2001a) for the intensity maps, reflect- 
ing the density structure. The A-variance spectrum in the 
intensity maps gradually steepens from a slope £ — 2.6 
on large scales to £ = 3.3 on the smallest scales. However, 
this behaviour is only visible when combining results from 
multiple instruments with different resolutions. The maps 
taken with each individual instrument appear to show only 
a power law within their limited dynamic range. The size- 
linewidth relation and structure function show more con- 
stant slope, as is expected from their lesser sensitivity to 
the driving scale. 

The velocity centroid PDFs show nearly Gaussian 
wings, as do the average line profiles, though with differ- 
ent widths. The PDFs of velocity difference as a function 
of lag, on the other hand, show non-Gaussian behaviour, 
with kurtosis values in some cases exceeding the exponen- 
tial value of six at the smallest scales, indicating correlated 
motion on these scales. 



We discuss below the physical implications of these 
observations, as deduced from our models. 



6.2. Model properties 

The clearest result we can draw from the models is that 
self-similar, power-law behaviour in both velocity and den- 
sity structure can be found between the dissipation scale 
and the driving scale even for the highly compressible tur- 
bulence that we model here. The behaviour is less clear 
but still definitely present after projection into 2D. The 
strength of the magnetic field also does not strongly in- 
fluence this conclusion, although stronger magnetic fields 
modify the slope of the density spectrum and do intro- 
duce measurable anisotropy into the turbulence, as shown 
in Fig. fH 

Above the driving scale, a flattening or turnover of the 
spectrum is apparent in all the measures we have studied, 
though it is most pronounced in the A-variance spectra. 
Conversely, numerical dissipation in the models causes a 
clear steepening of the spectrum as structure at smaller 
scales disappears. Although physical dissipation will not 
have the same detailed properties as numerical dissipa- 
tion, the general steepening of the spectrum will certainly 
occur. Similar behaviour is observed in incompressible tur- 
bulence (Lesieur 1997). 

The velocity centroid PDFs have kurtosis values in the 
range 2-4. The lowest values are distinctly sub-Gaussian, 
and appear to be produced onl y by driving at the largest 
available scales. We agree with Klessen (2000) that these 
non-Gaussian PDFs may be due to undersampling of a 
Gaussian distribution. No other physics that we have in- 
troduced produces such deviations. Observations of sub- 
Gaussian velocity centroid PDFs thus is consistent with 
driving from scales larger than the observed region. 

The models driven at the largest scale show marked 
super-Gaussian behaviour in their velocity difference 
PDFs at short lags, as reflected in kurtosis values exceed- 
ing even the exponential value of six. The occurrence of 
such values is thus a strong indication that the primary 
driving occurs on scales as large as or larger than those in 
the map. 



6.3. Turbulence properties 

Comparison between the observations of the Polaris Flare 
and the turbulence simulations constrains the mecha- 
nisms driving the turbulence in this cloud. Any mecha- 
nism that drives at an intermediate length scale, such as 
jets from embedded protostars, should produce character- 
istic features in the A-variance at that scale, which are 
not observed. The approximately self-similar, power-law 
behaviour seen in the observations is best reproduced by 
models where the energy is injected at large scales and 
dissipated at the small scales. The slight steepening of the 
A-variance spectra of the intensity structure, which seems 
also to be present in the velocity structure, is consistent 
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with turbulence models where the turbulence was recently 
injected and is decaying now due to the dissipation (see 
paper I). 

The large driving scale argues against protostellar out- 
flows being the main driving mechanism. An attractive 
alternative is driving by interactions with superbubbles 
and field supernova remnants (e.g. Mac Low et al. 2001, 
Avillez et al. 2000, Norman & Ferrara 1996). The Polaris 
Flare molecular cloud lies in the wall of a large cylindrical 
structure representing one of the nearest H I supershells, 
the North Celestial Polar Loop (Meyerdierks et al. 1991), 
adding additional support to this proposal. 

The kurtosis of the velocity difference PDFs in the 
Polaris Flare observations, especially at small and inter- 
mediate scales, reaches values close to 10, arguing for driv- 
ing from scales larger than the maps on which the high 
kurtosis appears, adding further support to the interpre- 
tation given above. We note that many of the actively star- 
forming regions observed by Miesch et al. (1999) have even 
higher kurtosis values at small scales, suggesting that the 
additional physics of strong self-gravity and local heating 
may produce additional effects that must be examined in 
future work. 

The dominant physical mechanism for dissi- 
pation in molecular cloud s was first shown by 



Zweibel fc Josafatsson (1983) to be ambipolar diffu- 



sion. Klcssen et al. (2000) showed that the length scale 
on which ambipolar diffusion will become important can 
be found by examining the ambipolar diffusion Reynolds 
number 



Ra = M A Liy n i/vA 



defined by 


Balsara (1996) 


Zweibel & Brandenburg (1997)|, where L 



(7) 



and 



ions, and v\ 



are the characteristic length and Alfven Mach number, 
jpi is the rate at which each neutral is hit by 
B 2 /Anp approximates the effective Alfven 
speed in a mostly neutral region with total mass density 
p = Pi + p n and magnetic field strength B. The coupling 
constant depends on the cross-section for ion-neutral 
interaction, and for typical molecular cloud conditions 
has a value of 7 « 9.2 x 10 13 cm 3 s _1 g _1 (e.g. Smith & 
Mac Low 1997). 

Setting the ambipolar diffusion Reynolds number 
Ra = 1 yields a diffusion length scale of 



Ld = V A /M.AVni 

w (0.041pc)A4 j4 



B 



10 pG 



10- 



10 J cm 



(8) 

S/2 

(9) 



with the ionization fraction x = pi/p n and the neutral 
number density n n — p n / P, with p = 2.36to#. If the ion- 
ization level in the Polaris Flare is low enough and the field 
is high enough, this length scale of order 0.05 pc would be 
directly resolved in the IRAM observations. We cannot yet 
unambiguously say whether the observed velocity spectra 
show a steepening there, similar to the downturn at the 



dissipation scale in the numerical models. If better obser- 
vations do show such a downturn in the future, that will 
be an indication of the dissipation scale. 

7. Summary 

We have applied several methods to characterise the ve- 
locity structure observed in the Polaris Flare molecular 
cloud over scales ranging from 0.015 pc to about 20 pc. 
We then applied the same methods to a large suite of 
computational hydrodynamic and MHD models of super- 
sonic, isothermal turbulence. The comparison between the 
observations and models with different parameters allows 
us to draw conclusions both about the properties of the 
analysis methods and the physical state of the observed 
region: 

— By measuring the average variation of velocity cen- 
troids as a function of the size of a virtual scanning 
beam in the observations and in simulated observa- 
tions, we recover Larson's (1981) size-linewidth rela- 
tion at scales where the turbulence shows an inertial 
range. This provides a method for measuring this rela- 
tion that does not rely on the identification of isolated 
clumps of gas. Comparing this variation with the cor- 
responding variation of the average total linewidth in 
the virtual beams allows estimation of the depth of the 
cloud along the line of sight. 

— The spectrum of the A-variance, a multi-dimensional 
wavelet transform ( |Stutzki et al. 199S ), shows devia- 
tions from inertial scaling behaviour at the scales of 
driving, dissipation, observational noise, and the tele- 
scope beam more clearly than other methods that we 
tested, because its effective filter function is better con- 
fined in the spatial frequency domain. However the 
method currently lacks intensity weighting, so that it is 
not reliable on observed maps with large empty areas. 

— The structure function characterizing the second mo- 
ments of the velocity difference PDFs duplicates infor- 
mation recovered by the scanning-beam size-linewidth 
relation or the A-variance of centroid velocities. The 
fourth moments measure correlated motion within a 
map. 

— We compare two proposed methods for measuring ve- 
locity PDFs by comparing their results to the known 
PDFs in our models. The average line profile, deter- 
mined in an optically thin line, is a better measure than 
the distribution of line centroid velocities. The width 
ratio of the two distributions depends on the depth of 
the observed region, as well as the Mach number of the 
flow, allowing us to infer a Mach number exceeding 10 
for the Polaris Flare. Sub-Gaussian velocity centroid 
PDFs are mainly produced by the cosmic variance in 
the case of large-scale driving. The influence of mag- 
netic fields and driving strength appears weak. 

— The power spectrum slope £ can be derived from the 
slope of the A-variance spectrum obtained for the ve- 
locity centroid maps of the Polaris Flare. At inter- 
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mediate scales in the KOSMA observations we find 
a slope in between 0.36 and 0.46, corresponding to 
C = 3.7. ..3.9, while at small scales in the IRAM map 
we find C = 3.6 . . . 4.4. Models of driven, supersonic 
turbulence show inertial range slopes of C =3.9-4.2, 
with slopes steepening to £ > 5 in the dissipation 
range. (Dissipation is due to numerical viscosity in 
these models.) The large error bars in the observed 
velocity spectra do not yet allow to draw a definite 
conclusion on a steepening in the dissipation range as 
seen in the models nor on a gra dual steepening across 
the full range of scales as seen by Bensch ct al. (2001a 
in the intensity maps of the Polaris Flare. Future ob- 
servations could do this. 

— The observations show super-Gaussian velocity differ- 
ence PDFs at small scales. Only driving at the largest 
scales in our models produces strongly super-Gaussian 
velocity difference PDFs at small lags. 

— The observed structure is consistent with hydrody- 
namic or MHD supersonic turbulence showing a com- 
plete spectrum from a driving scale at larger than 
10 pc, through an inertial range, to a dissipation scale 
under 0.05 pc. The combination of the density and the 
velocity scaling behaviour indicates a medium domi- 
nated by shocks at intermediate length scales creating 
thin sheets or filaments. Ambipolar diffusion could ex- 
plain the dissipation. 

— The main uncertainties come from the noise in the ob- 
servational data leading to large error bars in the ve- 
locity scaling functions and the unknown anisotropy 
of the velocity field. Anisotropic velocity fields could 
be produced by large-scale driving by supernovae or 
by strong magnetic fields but cannot be detected in 
observations bound to the plane of the sky. 
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